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A  multi-stepping  state  transition  matrix  approach  for  solving  a  linear  two-point 
boundary  value  problem  is  developed.  The  algorithm  employs  partitioned  state  transition 
matrix  of  the  Hamiltonian  system,  and  is  computationally  less  expensive  than  backward 
integration  of  differential  Riccatti  equation.  This  fact  makes  it  ideally  suited  for  online 
implementation.  The  application  of  this  technique  is  illustrated  for  a  finite  interval  moving 
mass  actuated  missile  guidance-autopilot  for  target  interception.  A  combination  of  feedback 
linearization  and  the  multi-stepping  linear  boundary  value  solution  algorithm  is  employed  in 
the  example.  Closed  loop  simulation  results  are  given. 


I.  Introduction 

Two-point  boundary  value  problems  result  from  the  application  of  optimal  control  theory  to  missile  guidance 
problems.  A  special  two-point  boundary  value  problem  of  interest  results  from  a  linear  dynamic  system  while 
optimizing  a  quadratic  performance  index  consisting  of  states  and  controls1  ’ 3.  As  an  example,  these  problems  arise 
in  the  derivation  of  guidance-control  laws  that  minimize  the  terminal  miss  distance,  while  penalizing  the  control 
effort.  Different  techniques  for  solving  the  two-point  boundary  value  problem  are  discussed  in  Reference  1. 
Techniques  such  as  the  shooting  method  require  the  solution  to  the  initial  value  problem  using  either  numerical 
forward  integration  of  the  differential  equations  or  the  use  of  state-transition  matrix  solution.  The  numerical 
integration  of  differential  equations  is  a  time  consuming  method,  while  the  state -transition  matrix  approach  suffers 
from  numerical  difficulties1  for  large  time  intervals  and  ill-conditioned  Hamiltonian  matrices.  Control  computation 
for  a  finite-interval  LQR  problem  can  also  be  posed  as  a  solution  to  the  Riccatti  differential  equation.  However, 
integrating  the  differential  Riccatti  equation  at  each  instant  of  time  backwards  may  not  feasible  for  real-time  control 
computation.  Off-line  gain  computation  and  implementation  requires  large  memory  from  the  on-board  processor. 
Moreover,  the  solution  may  not  be  useful  in  a  dynamic  setting  where  the  boundary  conditions  keep  changing  with 
time. 

An  analytical  state  transition  matrix  based  solution  has  been  discussed  in  Reference  1.  However,  the  state- 
transition  matrix  can  be  difficult  to  compute  for  large  time  periods.  This  paper  addresses  the  numerical  ill 
conditioning  problem  by  dividing  the  time-interval  into  multiple  intervals  and  employing  the  transition  matrix 
solution  in  each  subinterval.  This  approach  dramatically  improves  the  numerical  condition  of  the  problem,  while 
avoiding  the  need  for  numerical  integration.  One  of  the  byproducts  of  this  algorithm  is  the  solution  to  the  differential 
Riccati  equation.  The  technique  can  be  implemented  using  linear-algebraic  operations  available  in  software 
packages  such  as  LAPACK4.  The  numerical  algorithm  is  described  in  the  Section  II.  Section  III  describes  the 
application  of  this  technique  to  the  development  of  a  finite-interval  guidance-control  system  for  a  kinetic  warhead. 
Engagement  simulation  results  for  a  moving  mass  actuated  kinetic  warhead  are  presented  in  section  IV.  Section  V 
summarizes  the  conclusions  from  the  present  research. 


*  Research  Scientist,  Optimal  Synthesis  Inc.,  868  San-Antonio  Road,  Palo  Alto,  CA,  94303 
1  Chief  Scientist,  Optimal  Synthesis  Inc.  Associate  Fellow,  AIAA. 

*  Research  Scientist,  Optimal  Synthesis  Inc.,  Senior  Member 

§  Senior  Guidance  and  Control  Engineer,  Code  G23.  Associate  Fellow,  AIAA. 

©  Copyright  2005  by  Optimal  Synthesis  Inc.  All  Rights  Reserved. 

1 

Published  by  the  American  Institute  of  Aeronautics  and  Astronautics  with  Permission 


Report  Documentation  Page 

Form  Approved 

OMB  No.  0704-0188 

Public  reporting  burden  for  the  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering  and 
maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of  information, 
including  suggestions  for  reducing  this  burden,  to  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington 

VA  22202-4302.  Respondents  should  be  aware  that  notwithstanding  any  other  provision  of  law,  no  person  shall  be  subject  to  a  penalty  for  failing  to  comply  with  a  collection  of  information  if  it 
does  not  display  a  currently  valid  OMB  control  number. 

1.  REPORT  DATE 

2005 

2.  REPORT  TYPE 

3.  DATES  COVERED 

4.  TITLE  AND  SUBTITLE 

5a.  CONTRACT  NUMBER 

Multi- Stepping  Solution  to  Linear  Two  Point  Boundary  Value  Problems 

iv/t: — 

5b.  GRANT  NUMBER 

111  iYHSMlC  lllicgl  aicu  ^UllllUl 

5c.  PROGRAM  ELEMENT  NUMBER 

6.  AUTHOR(S) 

5d.  PROIECT  NUMBER 

5e.  TASK  NUMBER 

5f.  WORK  UNIT  NUMBER 

7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 

Naval  Surface  Warfare  Center, Dahlgren,VA, 22448 

8.  PERFORMING  ORGANIZATION 

REPORT  NUMBER 

9.  SPONSORING/MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 

10.  SPONSOR/MONITOR'S  ACRONYM(S) 

11.  SPONSOR/MONITOR'S  REPORT 
NUMBER(S) 

12.  DISTRIBUTION/AVAILABILITY  STATEMENT 

Approved  for  public  release;  distribution  unlimited 

13.  SUPPLEMENTARY  NOTES 

The  original  document  contains  color  images. 

14.  ABSTRACT 

see  report 

15.  SUBIECT  TERMS 

16.  SECURITY  CLASSIFICATION  OF: 

17.  LIMITATION  OF 

ABSTRACT 

18.  NUMBER 

OF  PAGES 

14 

19a.  NAME  OF 

RESPONSIBLE  PERSON 

a.  REPORT 

unclassified 

b.  ABSTRACT 

unclassified 

c.  THIS  PAGE 

unclassified 

Standard  Form  298  (Rev.  8-98) 

Prescribed  by  ANSI  Std  Z39-18 


Presented  at  the  AIAA  Guidance,  Navigation  and  Control  Conference,  August  15  -18,  2005,  San  Francisco, 
California. 


II.  Multi-Stepping  Technique 

The  two-point  boundary  value  problem  resulting  from  control  computation  of  a  linear  dynamic  system  while 
optimizing  a  quadratic  performance  index  can  be  defined  as  follows: 


z  =  Az  +  Bu 

min  J  =  zTfS  fzf  +  J  ( zT  Qz  +  uTRu)dt 

o 

subject  to  z(  0)  =  z0  and  z(tf)  =  zf 


(1) 


The  optimal  control  u  for  the  above  problem  is  a  function  of  the  co-state  vector  at  the  initial  A  (0),  which  can  be 
obtained  by  solving  the  following  two-point  boundary  value  problem: 


z 

~A  -BR'Bt~ 
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_A\ 

to 

1 

1 _ 

A 

,  z(0)  =  z0,  A(t  f)  =  S fz(t f) 


(2) 


The  coefficient  matrix  in  the  above  linear  system  of  differential  equations  is  known  as  the  Hamiltonian  matrix. 
The  solution  to  the  above  two-point  boundary  value  problem  can  be  obtained  if  the  state-transition  matrix  of  the 
Hamiltonian  matrix  can  be  computed  for  time  ( t  =  tj).  However,  in  most  problems,  the  state  transition  matrix  can  be 
difficult  to  compute  for  large  values  of  t  if  the  Hamiltonian  matrix  has  eigen-values  with  positive  real  parts.  The 
state  transition  matrix  can  also  be  difficult  to  compute  for  small  values  of  t  for  an  aggressive  choice  of  Q  and  R, 
which  can  result  in  an  ill-conditioned  Hamiltonian  matrix.  However,  for  sufficiently  small  t  the  state  transition 
matrix  can  be  computed  for  any  choice  of  Q  and  R.  This  fact  forms  the  basis  of  the  multiple  stepping  algorithm 
developed  in  this  paper. 

The  interval  [ 0  tj\  is  divided  into  sub-intervals  of  equal  length  ‘//’  such  as  [ 0  tj\,  [ti  2t/\, . |  (n-l)ti  /, | .  The 

length  of  the  interval  is  chosen  such  that  the  state  transition  matrix  be  computable  at  t  —  t/.  The  state  transition 
matrix  for  each  of  these  intervals  is  represented  as  shown  below: 

jY  =  eMt‘  =  eM(,2~h)  =  eM{h+t'~t')  =  eMt‘  = .  (3) 

where  M  is  the  Hamiltonian  matrix,  ti  and  t?  are  the  time  instants  at  the  end  of  the  first  and  second  intervals 
respectively.  The  solution  for  state  and  co-state  vector  at  end  of  each  time  instant  can  be  represented  as  shown 
below: 


zi 

zo 

~N 

0 

0 

0 

0" 

A0 

Z2 

= 

0 

N 

0 

0 

0 

z> 

0 

0 

0 

0 

N 

Z/ 

Z/-1 

aJ 

A-i 

The  central  idea  of  the  multi-stepping  approach  is  to  solve  for  the  initial  condition  on  the  co-state  vector 
recursively  using  this  matrix  equation.  The  multi-stepping  strategy  for  solving  the  initial  condition  on  the  co-state 
vector  is  demonstrated  in  a  two  interval  setting  below: 
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The  last  two  equations  ( A2  and  Z2  )  of  the  above  system  of  equations  can  be  expanded  as  shown  below: 

z2  =  Nnzx  +  NX2AX 
A2  =  N2xzx  +  N  22  A 1 


(6) 


where  Nfs  are  sub-matrices  of  the  state  transition  matrix  N.  Invoking  the  boundary  condition  A2  =  S  rZ2,  the 
following  equations  are  obtained: 

A  =[AT22 -SfNnWSfNn  -  Nlx\zx  =  Sxzx  (7) 

Si  in  the  above  expression  can  be  evaluated  since  the  state  transition  matrix  N  can  be  evaluated  and  therefore,  its 
sub-matrices.  Using  the  above  result  and  solving  the  system  of  equations  for  Ax ,  Zj  in  the  same  fashion  the 
following  expression  for  A0  is  obtained 

A0=[N22  -S1N12r1[S1Nn-N21]z0=S0z0  (8) 

The  above  expression  offers  a  computable  solution  to  the  initial  condition  on  the  co-state  vector.  The  multi 
stepping  strategy  obtains  the  relation  between  the  co-state  vector  A  and  the  state  vector  z  at  the  end  of  each  interval 

starting  with  Sf  at  t  =  tf.  It  is  well  known  that  the  matrix  connecting  these  two  vectors  is  the  solution  to  the 
differential  Riccatti  equation.  Therefore,  Sj  is  the  solution  of  the  finite  interval  Riccatti  equation  solution  at  time 
t=t].  As  the  number  of  time  intervals  is  increased  the  following  recursive  relation  can  be  used  to  compute  the 
Riccatti  equation  solution(S)  at  different  instances  of  time  in  between. 


A-i  =[^22 -N2l],  S(t  =  tf)  =  Sf 


(9) 


and  A0  =  S0Z0 


u(0)  =  -R-'BtA(0).  (10) 

The  multi-stepping  strategy  is  computationally  much  more  efficient  than  the  numerical  backward  integration  of 
differential  Riccatti  equation  at  each  instant  of  time,  since  it  uses  the  state  transition  matrix  to  propagate  states  and 
co-states.  It  also  offers  a  solution  when  conventional  techniques  requiring  solution  to  the  initial  value  problem  fail.  It 

should  be  noted  that  online  implementation  of  the  above  technique  only  requires  the  storage  of  the  matrix  N  on  the 
onboard  computer. 

In  the  following  two  sections  the  above  technique  will  be  applied  to  a  finite-duration  missile  guidance  problem. 

III.  Missile  Guidance 

Conventional  missile  control  is  typically  accomplished  using  actuators  like  fins  or  reaction  jets.  Moreover  the 
guidance  and  control  problems  are  decoupled  and  addressed  separately.  The  guidance  algorithm  generates  the 
acceleration  commands  for  the  autopilot.  The  autopilot  then  tracks  the  guidance  commands  to  achieve  target 
interception.  Recent  research5  has  advanced  a  technique  for  deriving  integrated  guidance-autopilot  system  in  a  finite 
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interval  setting.  This  formidation  of  the  integrated  guidance-autopilot  problem  forms  the  basis  for  the  present 
example. 

A  novel  moving  mass  based  actuation  system  was  proposed  in  Reference  6.  These  actuators  are  completely 
enclosed  within  the  envelope  of  the  kinetic  warhead  (KW).  An  integrated  approach  to  the  guidance  and  control  of 
the  KW  was  also  demonstrated  in  that  paper.  A  combination  of  continuous  time  feedback  linearization  and  pole 
placement  technique  was  used  for  the  integrated  guidance-control  system  design  (IGCSD).  The  controller’s  task  was 
to  regulate  the  instantaneous  line  of  sight  rate  of  the  target  with  respect  to  the  missile.  A  commercially  available 
nonlinear  control  system  design  software7'8  was  been  used  for  controller  design. 

An  IGCSD  formulation  that  is  suitable  for  finite  duration  implementation  will  be  developed  in  this  section. 
Figure  1  shows  a  schematic  of  the  kinetic  warhead’s  target  interception  scenario.  Interception  occurs  when  the 
relative  position  vector  between  the  KW  and  the  target  goes  to  zero.  Therefore,  the  control  objective  is  to  drive  the 
terminal  relative  position  vector  between  the  KW  and  target  to  zero.  However,  only  two  components  of  the  relative 
vector  are  controllable.  The  component  of  the  relative  position  vector  along  the  longitudinal  body  axis  is  not 
controllable  using  the  moving  mass  actuators.  This  leaves  the  mass  movement  along  the  body  y-axis  and  the  mass 
movement  along  the  body  z-axis  as  the  control  variables  in  the  problem. 


Target 


Figure  1  :  Kinetic  Warhead  and  Target  Engagement  Scenario 

The  control  is  achieved  by  generating  acceleration  components  normal  to  the  velocity  vector.  The  interception 
problem  is  posed  as  a  two  parameter  control  problem  in  a  plane  normal  to  the  line  of  sight  vector.  In  this  example 
the  relative  vector  components  along  the  y  and  z  directions  of  a  LOS  frame  will  be  controlled  over  a  finite  duration. 
The  LOS  frame  is  chosen  such  the  x-axis  coincides  with  the  line  of  sight  vector  at  the  initial  time.  The  components 
of  the  initial  relative  position  vector  A  v  and  A .  are  zero  in  this 
and  z  actuators  are  shown  below: 

5yc  — >  Uy  — >  5y  ^  ^  ^  ^  ^  ^  A  ^ 

3-c  uz  q  w  - ►  A, 


frame  of  reference.  The  control  chains4  for  the  y 


Av 

A_ 


(11) 
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8  yc  -  Position  Command  to  the  y  actuator 
8ZC  -  Position  Command  to  the  z  actuator 

Uy:  =  kp  (Svc  —S')  —  kvSv  +  m  av  -  Force  applied  on  the  y  actuator, 
kp  and  kv  are  the  position  servo  gains 

U,  =  kp  (8  zc  —  8 z)  —  kv8,  +  rn_az  -  Force  applied  on  the  z  actuator 

8  -  Actual  Position  of  the  y  actuator 

y 

8 _  -  Actual  Position  of  the  z  actuator 

ay  -  Acceleration  component  along  the  body-frame  y-axis 
a~  -  Acceleration  component  along  the  body-frame  z-axis 

my  -  mass  of  the  y-actuator 
mz  -  mass  of  the  z-actuator 

q  -  Pitch  rate  of  the  KW 
r  —  Yaw  rate  of  the  KW 

v  -  Body-frame  velocity  component  along  the  y-axis 
w  -  Body-frame  velocity  component  along  the  z-axis 

The  variables  A  v ,  A z  are  the  components  of  the  relative  velocity  vector  in  the  inertial  frame. 


The  equations  of  motion  of  the  KW  given  in  Reference  6  are  used  for  simulation  and  controller  design.  The 
simulation  involves  a  six-degree  of  freedom  KW  dynamic  model  and  two  two-degree  of  freedom  actuator  models 
for  the  y  and  z  moving  masses.  A  three  degree  of  freedom  particle  model  is  used  to  simulate  the  dynamics  of  the 
target.  North-East-Down  (NED)  inertial  frame  of  reference  is  used  to  represent  the  position  of  the  KW  and  the 
target.  Feedback  linearization  is  performed  at  each  instant  of  time  numerically  using  the  nonlinear  control  system 
design  software.  The  resulting  dynamical  system  consists  of  two  sixth-order  integrators  along  each  control  chain 
with  a  pseudo-control  inputs  upy  and  up:\ 

A<6)  =  U  &  A<6)  =  U  (12) 

y  py  z  pz  V  ' 

The  target  interception  problem  can  now  be  posed  as  a  finite-duration  linear  quadratic  control  problem  on  these 
pure-integrator  dynamical  systems,  with  a  quadratic  performance  index  consisting  of  the  final  values  of  the  relative 
position  components  and  the  integral  of  a  quadratic  form  in  the  pseudo-control  variables. 

The  resulting  optimal  control  problem  can  be  solved  using  the  multi-stepping  algorithm  discussed  in  Section  2. 
The  solution  will  be  in  terms  of  time-to-go.  An  approximate  value  of  the  time-to-go  can  be  computed  using  the 
instantaneous  values  of  the  relative  velocity  and  position  vectors  between  the  target  and  the  kinetic  warhead.  Figure 
2  shows  a  block  diagram  of  the  closed-loop  system.  Engagement  results  using  the  finite-duration  controller  are  given 
in  the  following  section. 
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NST 


Position  Force 


Figure  2.  Block  Diagram  of  the  Controller  Implementation 


IV.  Engagement  Simulations 

The  following  design  parameters  were  chosen  for  the  controller: 

Q  =  100*diag([0  01001  00]) 

R  =  1 

Sf=  zeros(6,  6);  S/1,1)  =  le6 
n_steps  =  20 

The  parameter  n_steps  =  20  refers  to  the  number  of  intervals  employed  by  the  multiple  stepping  approach  to 
compute  the  finite  interval  control. 

A.  Engagement  Scenario  #1 

Initial  conditions  on  the  attitude  of  the  KW  are  chosen  to  result  in  zero  angle  of  attack  and  zero  angle  of  side-slip 
at  time  t  =  0.  This  can  be  done  by  choosing  the  pitch  attitude  angle  same  as  the  flight  path  angle  and  the  yaw  attitude 
angle  same  as  the  heading  angle  of  the  KW.  The  roll  attitude  angle  has  always  been  chosen  as  zero.  The  body  rates 
along  all  three  axes  have  also  been  chosen  as  zero  at  time  t  =  0.  Initial  displacements  of  the  moving  masses  and  their 
initial  speeds  have  also  been  set  to  zero. 

Shown  in  Figure  3,  Figure  4,  and  Figure  5  are  the  trajectories  of  the  KW  and  the  target  in  3D,  vertical  plane  and 
horizontal  plane  respectively.  The  target  in  this  case  is  descending  from  a  higher  altitude  and  is  initially  located 
south-east  of  the  KW.  The  finite  duration  IGCSD  minimizes  the  y  and  z  components  of  the  relative  position  vector 
at  the  final  time  as  shown  in  Figure  6.  The  initial  values  of  the  y  and  z  components  are  zero  in  the  inertial  frame 
owing  to  the  definition  of  the  frame.  Flowever,  they  keep  changing  due  to  the  difference  in  the  velocity  vector  and 
the  LOS  vector  before  assuming  a  very  small  terminal  value.  A  terminal  miss-distance  of  0.26481  (ft)  which  is  less 
than  the  diameter  of  the  KW  resulted  indicating  successful  target  interception  by  the  KW.  The  actual  and 
commanded  displacements  of  the  y  and  z  actuators  are  shown  in  Figure  7  and  Figure  8.  The  commanded 
displacements  are  the  result  of  feedback  linearization  and  finite  duration  control.  The  y  moving  mass  traverses  most 
of  its  stroke  length  whereas  the  z  moving  mass  travels  smaller  distance.  The  magnitude  of  the  mass  displacement  is 
governed  by  the  acceleration  requirements  along  the  y  and  z  directions.  After  the  first  few  seconds  of  transient  the 
attitude  of  the  KW  remains  constant  and  the  body-rates  go  to  zero  as  shown  in  Figure  9.  It  should  be  noted  that  all 
plots  are  generated  using  non-dimensional  variables. 
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3-D  Trajectories  of  the  KW  and  the  Target 

- Warhead 

- Target 


Figure  3.  3D  Interception  Trajectories 


Vertical  Plane  Trajectories  of  the  KW  and  the  Target 


Figure  4.  Vertical  Plane  Trajectories 
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Figure  5.  Horizontal  Plane  Trajectories 


Y  &  Z  Components  of  the  Relative  Position  Vector 


Figure  6.  Components  of  the  Relative  Position  Vector  in  LOS  Frame 
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Position  of  the  Y-Actuator  Mass 


Figure  7.  Displacement  of  the  Y-Actuator  Moving  Mass 


Position  of  the  Z-Actuator  Mass 


Figure  8.  Displacement  of  the  Z-Actuator  Moving  Mass 


9 

Published  by  the  American  Institute  of  Aeronautics  and  Astronautics  with  Permission 


Presented  at  the  AIAA  Guidance,  Navigation  and  Control  Conference,  August  15  -18,  2005,  San  Francisco, 
California. 


KW  Body-Rates 


Figure  9.  Body  Rates  of  the  KW 


B.  Engagement  Scenario  #2 

The  target  in  this  scenario  is  approaching  from  the  south-west  direction  of  the  initial  position  of  the  KW.  NED 
frame  trajectories  of  the  KW  and  the  target  are  shown  in  Figure  10,  Figure  11,  and  Figure  12.  Terminal  miss- 
distance  obtained  from  the  simulation  is  0.22489(ft)  indicating  successful  interception. 

3-D  Trajectories  of  the  KW  and  the  Target 

- Warhead 

- Target 


Figure  10.  3D  Interception  Trajectories 
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Vertical  Plane  Trajectories  of  the  KW  and  the  Target 


North  Position 


Figure  11.  Vertical  Plane  Trajectories 


Horizontal  Plane  Trajectories  of  the  KW  and  the  Target 

0.1  - 1 - 1 - 1 - 1 - 1 -  "i  '  - 

- Warhead 

- Target 


0.05 


-0.25 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 

-0.12  -0.1  -0.08  -0.06  -0.04  -0.02  0  0.02 


East  Position 


Figure  12.  Horizontal  Plane  Trajectories 


Components  of  the  relative  vector  in  the  LOS  frame  are  shown  in  Figure  13  which  start  and  end  at  zero. 
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Y  &  Z  Components  of  the  Relative  Position  Vector 


Figure  13.  Components  of  the  Relative  Position  Vector  in  the  LOS  Frame 


Position  of  the  Y-Actuator  Mass 


Figure  14.  Displacement  of  the  Y-Moving  Mass 


12 

Published  by  the  American  Institute  of  Aeronautics  and  Astronautics  with  Permission 


Presented  at  the  AIAA  Guidance,  Navigation  and  Control  Conference,  August  15  -18,  2005,  San  Francisco, 
California. 


Position  of  the  Z-Actuator  Mass 


Figure  15.  Displacement  of  the  Z-Moving  Mass 


V.  Conclusion 

A  multi-stepping  approach  was  developed  to  solve  linear  two-point  boundary  value  problems  arising  missile 
guidance-control  problems.  This  approach  is  computationally  less  expensive  than  the  backward  integration  of  the 
differential  Riccatti  equation.  The  algorithm  has  been  applied  to  a  finite-duration  integrated  guidance-control  system 
design  of  a  moving  mass  actuated  kinetic  warhead.  Feedback  linearization  is  used  to  convert  the  nonlinear  control 
problem  into  a  linear  control  problem.  Finite-interval  linear  optimal  control  techniques  are  then  used  to  pose  the 
control  computation  problem  as  a  solution  to  a  two-point  boundary  value  problem.  Multi-stepping  approach 
presented  in  this  paper  has  been  used  to  compute  the  control. 
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